/* 
   xmunipack - value units (intensity, flux, etc..) conversions
   
   Copyright © 1997-2011 F.Hroch (hroch@physics.muni.cz)
   
   This file is part of Munipack.
   
   Munipack is free software: you can redistribute it and/or modify
   it under the terms of the GNU General Public License as published by
   the Free Software Foundation, either version 3 of the License, or
   (at your option) any later version.
   
   Munipack is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.
   
   You should have received a copy of the GNU General Public License
   along with Munipack.  If not, see <http://www.gnu.org/licenses/>.

*/

/*

   This routine represents physical value intensity and provides
   its represenattion in various physical units.

   Currently only HST magnitude - photometry convention is used.

*/


#include "fits.h"
#include <wx/wx.h>

#ifdef __WXDEBUG__
#include <wx/debug.h>
#endif


using namespace std;

FitsValue::FitsValue(): type(UNIT_COUNT),hascal(false) {}

FitsValue::FitsValue(const FitsArray& a): type(UNIT_COUNT),hascal(false)
{
  arrays.push_back(a);
  hascal = init(arrays,photflam,exptime);
}

FitsValue::FitsValue(const FitsArray& r,const FitsArray& g,const FitsArray& b):
  type(UNIT_COUNT),hascal(false)
{
  arrays.push_back(r);
  arrays.push_back(g);
  arrays.push_back(b);
  hascal = init(arrays,photflam,exptime);
}

bool FitsValue::init(const std::vector<FitsArray>& arrays, 
	  std::vector<double>& photoflam, std::vector<double>& exptime)
{
  bool hc = true;

  for(vector<FitsArray>::const_iterator k=arrays.begin(); k!=arrays.end();++k){
    FitsHdu a(*k);

    wxString u,v;
    double ph,t;
    
    u = a.GetKey("PHOTFLAM");
    v = a.GetKey("EXPTIME");

    hc = hc && ( !u.IsEmpty() &&  u.ToDouble(&ph)) &&
               ( !v.IsEmpty() &&  v.ToDouble(&t));

    photflam.push_back(ph);
    exptime.push_back(t);
  }
  return hc;
}


void FitsValue::SetType(units_type t)
{
  type = t;
}

units_type FitsValue::GetType() const
{
  return type;
}

std::vector<double> FitsValue::Get(int i, int j) const
{ 
  std::vector<double> vals;

  for(vector<FitsArray>::const_iterator k=arrays.begin(); k!=arrays.end();++k){
    FitsArray a(*k);

    if( 0 <= i && i < a.Width() && 0 <= j && j < a.Height() )
      vals.push_back(a.Pixel(i,j));
    else
      vals.push_back(0.0);
  }
  return vals;
}

wxString FitsValue::Get_str(int i, int j) const
{
  wxString line;

  for(size_t k = 0; k < arrays.size(); k++) {
    FitsArray d = arrays[k];

    wxString a;
    if( 0 <= i && i < d.Width() && 0 <= j && j < d.Height() ) {
      
      if( type == UNIT_COUNT ) {
	if( d.Bitpix() > 0 )
	  a.Printf("%d",int(d.Pixel(i,j)+0.5));
	else
	  a.Printf("%.5g",d.Pixel(i,j));
      }
      else if( hascal ) {

	double z = photflam[k];
	double t = exptime[k];

	if( type == UNIT_MAG ) {

	  double f = d.Pixel(i,j)/t + z;
	  a.Printf("%.3f",f > 0.0 ? -2.5*log10(f) : 99.999);
	}
      }
    }
    if( ! a.IsEmpty() )
      line += " " + a;
  }
  return line;
}

wxString FitsValue::Label_str(int n)
{
  switch(n){
  case UNIT_COUNT:  return "counts (ADU)";
  case UNIT_MAG:    return "mag/arcsec²";
  case UNIT_PHOTON: return "photons/s/cm²";
  case UNIT_ERG:    return "erg/s/cm²";
  default:          return wxEmptyString;
  }
}

wxArrayString FitsValue::Label_str()
{
  wxArrayString a;

  for(int i = UNIT_FIRST+1; i < UNIT_LAST; i++)
    a.Add(Label_str(i));
	  
  return a;
}
